Hexatic and mesoscopic phases in the 2D quantum Coulomb system 
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We study the Wigner crystal melting in a two dimensional quantum system of particles interacting 
via the 1/r Coulomb potential. We use quantum Monte Carlo methods to calculate its phase 
diagram, locate the Wigner crystal region, and analyze its instabilities towards the liquid phase. 
We discuss the role of quantum effects in the critical behavior of the system, and compare our 
numerical results with the classical theory of melting, and the microemulsion theory of frustrated 
Coulomb systems. We find a Pomeranchuk effect much larger then in solid helium. In addition, 
we find that the exponent for the algebraic decay of the hexatic phase differs significantly from the 
Kosterilitz-Thouless theory of melting. We search for the existence of mesoscopic phases and find 
evidence of metastable bubbles but no mesoscopic phase that is stable in equilibrium. 
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PACS numbers: 29.25.Bx. 41.75.-i, 41.75.Lx 

The Wigner crystal (WC) melting has been a subject 
of intense study over the years. [3, 0, Q A better under- 
standing of this process is particularly important in the 
two dimensional (2D) one component plasma (OCP) with 
1/r Coulomb interactions, since it could explain many 
features in systems such as electrons at interfaces, dusty 
plasmas, MOSFETs, and charged colloids. 

Although in 2D, a true crystalline order is not possi- 
ble at any finite temperature T |, a quasi- long range 
translational order is stable at strong coupling. Upon 
melting from the WC a variety of mechanisms have been 
3roposed.[l| At high temperature, Halperin and Nelson 
H suggested that the 2D melting is a two-step pro- 
cess, with a hexatic phase in between. Hexatic phases 
have been seen in colloids, liquid crystals, and dusty 
plasmas, @ and found in some classical simulations. @, [I| 
However, there is disagreement about its location, critical 
exponents, and order. The influence of quantum effects 
is unclear. 

In the ground state, Jamei et al. 0] argue by a mean- 
field approach that a first order WC-liquid transition at 
zero temperature is disallowed; instead the transition is 
mediated by microemulsion phases (e.g. stripes or bub- 
bles). This mechanism might be related to the stripe 
order in some strongly correlated materials and anoma- 
lous effects in experiments on MOSFETs. fiol] 

In this Letter, we examine a 2D system of N quan- 
tum distinguishable (Boltzmannon) particles that inter- 
act with a Hamiltonian: 
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where is the location of the i'th particle. The system 
is parameterized by two dimensionless parameters which 
we choose to be r s and T. We use a — r s aa as the unit 
of length, with ao the Bohr radius, Rydberg as the unit 
of energy, and temperature; the density is l/7r(r s ao) 2 . 



The system is neutralized by a rigid homogeneous back- 
ground, which gives rise to the constant Vt> g . Our goal is 
to map out the phase diagram, to locate the WC phase 
boundaries, and analyze the transition to the liquid. To 
do so, we use quantum Monte Carlo (QMC) techniques, 
which arc uniquely suited for calculations of the proper- 
ties of strongly coupled Coulomb systems. 

At finite temperature, we use path integral Monte 
Carlo (PIMC).[ll| to sample configurations from the 
thermal density matrix. In contrast to fermion systems, 
for a system of distinguishable particles, there is no sign 
problem. All systematic errors, such as the time step and 



finite size errors, can be studied with extrapolation. 22 1 
Our calculations involve approximately 100 phase points 
at system ranging from 200 to 2248 particles and between 
10 and 8000 slices in imaginary time. 

At T = 0, we use diffusion Monte Carlo (DMC), a pro- 
jector method that takes as input a trial wave function 
and filters from it, all excited states. Calculating 
DMC energies in this system is formally exact as long as 
has a non-zero overlap with the ground state. 

General properties of the phase diagram - Our pro- 
posed phase diagram is shown in Fig. [T] The boundary 
between the WC and the liquid phase presents a clear re- 
entrant behavior at low T. Although such behavior has 
been seen in helium, 14, for the 2D-OCP the effect 
is significantly stronger. We find the ratio between the 
minimal and maximal density of the crystal 0.6 compared 
to 0.03 for He 3 and 3- 10~ 5 for He 4 . The re-entrant liquid 
terminates in a zero temperature quantum critical point 
(QCP) at r s ~ 66.5 ± 0.2. This differs by 10 % from pre- 
vious DMC simulations; [12| our calculations did a more 
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careful extrapolation to the thermodynamic limit. 

At large r s we compare our results with simulations 
done for the classical OCP; classical simulations found a 
WC phase stable for r > 140 and a liquid phase stable 
for T < 120 ,8] where T = 2/r s T. Seen as lines in Fig.HJ 
these values compare favorably with our results at large 
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FIG. 1: Phase diagram of the quantum 2D Coulomb sys- 
tem calculated using QMC as obtained by analyzing corre- 
lation functions. Points indicate location of QMC calcula- 
tions. The phase boundaries are interpolations between these 
points. The dotted blue line indicates a rough separation be- 
tween the classical and the quantum regions: when the ther- 
mal de Broglie wave length equals the interparticle spacing. 
The solid straight blue lines are the classical results of Ref. [S| . 
The inset is a density profile of the metastable bubbles seen 
in PIMC. 



r s . For 120 < T < 140, we found a hexatic quasi-long 
range order; this order extends into the quantum regime 
by bending with respect to the classical line as the density 
is increased (i.e. for r s < 200). 

At r s ~ 60 the melting temperature decreases until 
it reaches a tricritical point (TCP), located at r s ~ 50 
and T ~ 80 fiRyd, where the isotropic liquid, the hex- 
atic phase, and the WC coexist. We have used the 
Clausius-Clapeyron relation to understand this behavior. 
For Coulomb systems it readsfl6j|: 



dT _ 2 SK + SU 
d{l/r s ) ~ F SU ' 



(2) 



where SU = U C — Uh and SK = K c — Kh are the changes 
in the internal and kinetic energies between the crystal 
(U c , K c ) and hexatic (Uh, Kh) phase. The slope is shown 
in Fig. [21 At the downturn, the slope is zero, while it 
diverges at the TCP, where the entropy difference, and 
hence SU, goes to zero. In fact, Eq. [2] shows that a first 
order line for any Coulomb system, must become second 
order at the nose of a re-entrant phase. 

There are two WC-to-liquid boundaries: the quantum 
melting line, going from the QCP to the TCP, and the 
thermal melting line, flowing from the high tempera- 
ture side of the TCP towards lower density (large r s ), 
and characterized by a two-step process mediated by a 
hexatic phase. The order of the phase transitions is es- 
tablished by examining the internal energy per particle 
U = U(T) at a series of densities r s e (48, 400) through 
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FIG. 2: dT/d(l/r 3 ) vs l/r s computed using the Clausius- 
Clapeyron relation, Eq. [2] for N = 562, compared with the 
results obtained directly from the phase diagram in Fig. [1] for 
the first order WC-to-hexatic boundary. 



a wide range of temperatures. The behavior at r s = 55 
is shown in Fig. [3J A second order (or weakly first or- 
der) phase transition is evident for the quantum melt- 
ing line, indicated by a kink in the internal energy when 
moving from the re-entrant liquid phase to the WC. In- 
creasing the temperature, U(T) shows a sharp jump at 
the boundary between the solid and hexatic phase, in- 
dicating a first order WC-hexatic transition. The inter- 
nal energy is continuous between the hexatic phase and 
the liquid, indicating a continuous hexatic-liquid transi- 
tion. This is in disagreement with the classical theory 
of melting, which suggests two continuous transitions of 
Kosterlitz-Thouless (KT) type. It should be noted that 
as one approaches the classical limit, r s (« 400), the 
jump in energy shrinks and becomes barely distinguish- 
able from a continuous transition. 
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FIG. 3: Internal energy per particle versus temperature at 
r s — 55 for N = 562. The inset plots both the number of 
defects and the internal energy scaled to arbitrary units and 
shifted. 



Thermal melting line - The thermal melting line is 
characterized by the presence of a hexatic intermediate 
phase with an order parameter: 
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where Vt'(r i ) = g 2m exp(6i0y); the summation is over 
nearest neighbor particles j surrounding i, [23} and Oij 
is the angle between an arbitrary fixed vector and the 
vector joining Ti to rj. The algebraic decay of ge(r) 
signals a quasi-long range orientational order and this, 
combined with a lack of translational order, defines the 
hexatic phase. Additionally, the orientational order of 
the hexatic phase manifests itself via a sixfold anisotropy 
in the structure factor S(k). Fig. |4] shows the behavior 
of the S(k) in the liquid, solid, and hexatic phases at 
r s = 55, and ge(r) as a function of temperature. There 
is a clear qualitative change in the decay of g§ (r) for 
156.2 iiRyd < T < 150.3 /J,Ryd, as we transition from 
the liquid to the hexatic phase. 

The KT theory predicts that for hexatic systems in 
the XY universality class g§(r) oc r -I?6 , with -q^ < 1/4, 
and 776 reaching its maximum 1 /4 at the transition to the 
liquid. Wc fit the large distance values of ge(r), to de- 
termine r]6 and found that, independent of temperature, 
it undergoes the hexatic-to- liquid transition for 776 — 2. 
The disagreement with the KT theory is likely due to 
the 1/r interaction. Indeed, critical exponents in other 
models, for example the Ising model, are known to be 
strongly affected by the range of the interaction. 
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Defects play an important role in 2D melting. We de- 
fine defects using a Voronoi construction of each particle's 
centroid. In a perfect hexagonal lattice, each centroid 
has six neighbors. Defective centroids, say with 5 or 7 
neighbors, are denoted as vortices and antivortices re- 
spectively. According to the KT theory the WC-hexatic 
transition is caused by the unbinding of dislocations (par- 
ticles defined as strongly bound vortex-antivortex pair) , 
while the hexatic-liquid transition is driven by the vortex- 
antivortex unbinding. However, important differences 
are found in our system. At the verge of the WC melting, 
instead of the expected continuous rise of dislocations, 
the number of defects jumps discontinuously. We see 
a strong correlation between the energy and the num- 
ber of defects; U(T) in the hexatic and solid phases is 
controlled entirely by their number. Therefore, the first 
order jump seen in U (T) at the WC-hexatic transition 
is due to the proliferation of dislocations in the system, 
which can be induced by grain boundaries, 18, 3| present 
in the WC phase. The inset of Fig. [3] shows this pro- 
portionality by plotting both the defect number and the 
energy of the various phases on the same axis, with the 
U(T) scaled to arbitrary units and shifted. We note that 
the first order transition in a long range system is uncon- 
ventional because the presence of the background does 
not allow macroscopic separation of phases with different 
densities. In the hexatic-liquid transition, the number of 
defects grows slowly while the energy of the system rises 
at a much faster rate suggesting the presence of an un- 
binding process possibly related to the vortex-antivortex 
melting. [5j] 

Quantum melting line - The quantum melting line 
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FIG. 4: Top: ge(r) at different temperatures for r 3 = 55, N = 
2248. The lower two temperatures show algebraic decay while 
the higher temperatures decay faster than algebraically. Mid- 
dle: rje (T) in the hexatic phase for N = 562 at various tem- 
peratures and densities (for r 3 values, see legend on figure). 
All exponents satisfy 776 < 2 (shown as dotted line). Bot- 
k y ) in a liquid, hexatic, and solid respectively for 



torn: S^kx, 

r s = 55, N 



562, temperatures in caption. 



goes from the QCP to the TCP. Jamei et al.\jj^ and Or- 
tix et al.^j^ showed recently by a mean-field Hamilto- 
nian that the system undergoes a WC-liquid transition 
at T — through a sequence of intermediate phases 
such as stripes, or other "microemulsion phases" . Their 
shape and size depend on the parameters of the mean- 
field Hamiltonian. At the critical density r°, where the 
free energy of the crystal equals that of the liquid, the 
optimal geometry is found to be alternating liquid and 
crystal stripes of width Wq — oexp(A), with a a lat- 
tice cutoff proportional to the mean interparticle spac- 
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47r 2 e 2 c7/A/^, where a the surface tension, 
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ing, A = 47r"e" 
A^i c = /icrys(^s) ~~ Miiq( r s) the chemical potential differ- 
ence between the crystal and liquid at criticality. 

At finite T, our PIMC results suggest a second order 
(or very weakly first order) direct transition. We care- 
fully searched for mesoscopic phases in the surrounding 
region. At r s rs 55, where the transition temperature to 
the re-entrant liquid is relatively high (T ~ 50 fiRyd), 
we find no evidence for mesoscopic phases. At larger 
r s and, hence, lower transition temperature, the PIMC 
simulations occasionally yield inhomogenous phases. A 
density profile of such a phase is shown in the inset of 
Fig. [TJ However, those inhomogenous structures appear 
to be metastable because their internal energy, measured 
in PIMC at low T (< 10 fiRyd) and DMC at T=0, is 
always higher than the homogenous phases. To clar- 
ify this situation, we performed additional DMC sim- 
ulations to determine the parameters of the mean-field 
model. In the liquid phase, we used a pair product or Jas- 
trow trial function #u q uid = ]li<j exp\-UR P A(\Tj - rj\)], 
with ujipa the RPA Jastrow function. [2l| For the solid 
and stripe phase we multiplied by a Gaussian, tieing 
the distinguishable particles to predefined sites, {1^} ar- 
ranged to correspond to a crystal or stripe: ^crystal = 
'I' liquid exp[- J2i a \ r i ~ h\ 2 ]- We fit the T = liquid and 
crystal energies to a polynomial and, assuming the tran- 
sition is weakly first order, obtain r c s = 66.5 ± 0.2 and 
A[i c = 65 ± 5 [iRyd. 

We estimate the surface tension a by computing the 
energy at r c s for different number of stripes, 24J and hence 
different surface length in the system: er s tripc = 1-55 ± 
0.09 nRyd/ciQ. Plugging cr s tripc and A/i c in the equation 
for A gives A « 3 x 10 5 , which leads to the width of the 
stripes at the critical point of Wq « 10 10 : much larger 
than any physical system. The large widths are a result 
of very small difference in the chemical potential between 
the homogeneous liquid and crystal phases at the QCP, 
characteristic of the 1/r interaction. Away from r c s , the 
periodicity of the alternating stripes will become even 
larger. Therefore, the instability towards "microemulsion 
phases" has a characteristic emergent length Wq so large 
that it cannot explain non-Fermi liquid behavior seen in 
2D experimental setups. (Iol| 

Conclusion - We have studied the WC melting in the 
2D quantum Coulomb system of particles with Boltz- 
mann statistics. Its phase diagram shows thermal and 
quantum melting. The thermal melting is mediated by a 
hexatic phase, but significant deviations from the classi- 
cal KT theory are found. The WC-hexatic transition is 
first order, driven by the proliferations of defects in the 
crystal, where they are assembled into grain boundaries. 
The hexatic parameter r]Q « 2 belongs to the universality 
class of long-range models and reveals the importance of 
the long-range interaction in determining the transition 
properties. In the low temperature region, we found a 
strong Pomeranchuck effect, much larger than for solid 



helium. At low temperatures, we did not find stable mi- 
croemulsion phases. An estimate of their size is exceed- 
ingly large, which makes this kind of phases impossible to 
see in any physical system if driven by correlation effects 
only. 
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